Consensus genomic regions associated with multiple abiotic stress tolerance in wheat and implications for wheat breeding

In wheat, a meta-analysis was performed using previously identified QTLs associated with drought stress (DS), heat stress (HS), salinity stress (SS), water-logging stress (WS), pre-harvest sprouting (PHS), and aluminium stress (AS) which predicted a total of 134 meta-QTLs (MQTLs) that involved at least 28 consistent and stable MQTLs conferring tolerance to five or all six abiotic stresses under study. Seventy-six MQTLs out of the 132 physically anchored MQTLs were also verified with genome-wide association studies. Around 43% of MQTLs had genetic and physical confidence intervals of less than 1 cM and 5 Mb, respectively. Consequently, 539 genes were identified in some selected MQTLs providing tolerance to 5 or all 6 abiotic stresses. Comparative analysis of genes underlying MQTLs with four RNA-seq based transcriptomic datasets unravelled a total of 189 differentially expressed genes which also included at least 11 most promising candidate genes common among different datasets. The promoter analysis showed that the promoters of these genes include many stress responsiveness cis-regulatory elements, such as ARE, MBS, TC-rich repeats, As-1 element, STRE, LTR, WRE3, and WUN-motif among others. Further, some MQTLs also overlapped with as many as 34 known abiotic stress tolerance genes. In addition, numerous ortho-MQTLs among the wheat, maize, and rice genomes were discovered. These findings could help with fine mapping and gene cloning, as well as marker-assisted breeding for multiple abiotic stress tolerances in wheat.

MQTLs associated with multiple abiotic stress tolerance. Only 2888 QTLs (out of 3102 total QTLs) were utilized for projection onto the consensus map, the remaining 214 QTLs could not be used for projection because of the lack of necessary information (required for projection) from the source studies. Only 2229 QTLs out of 2888 could be projected onto the consensus map, leaving 659 QTLs un-projected on the consensus chromosomes owing to either of the following reasons: (i) lack of common markers between the initial and consensus map, (ii) large CI of initial QTLs. Of 2229 projected QTLs, 2139 QTLs were clustered into 147 QTL hotspots, leaving 13 singletons and 47 unassigned to any hotspot. Out of 147 QTL hotspots, 13 hotspots included initial QTLs only from a single study, and therefore could not be considered as true MQTLs. The 134 true MQTLs were distributed on all 21 chromosomes, with a minimum of only one MQTL on chromosome 3A and a maximum of 10 MQTLs on chromosome 7A (Table S4;   The MQTLs had an average CI of 9.93 times smaller than the initial QTLs (the average CI of initial QTLs was 19.58 cM), and there were substantial differences in CI reduction among different chromosomes. On chromosomes 2D and 3B, the average CI of MQTL was reduced by 22.84 and 20.92 times, respectively, followed by 19.5 and 18.68 times on chromosomes 7A and 2A, respectively (Fig. 2c). Interestingly, only a single MQTL was identified on chromosome 3A, which showed a CI of only 0.01 cM, whereas, the initial QTLs located on chromosome 3A had an average CI of 16 MQTL validation and the availability of known genes within MQTL regions. There were considerable overlaps between the MQTLs discovered through meta-analysis and the significant markers (or markertrait associations; MTAs) identified through the GWAS approach that are associated with various abiotic stresses in the wheat genome. Remarkably, 817 MTAs identified in wheat GWAS for various abiotic stresses overlapped with 76 MQTLs out of the 132 physically anchored MQTLs (Table S5, Fig. 3). There were significant differences in the number of MTAs overlapping with an individual MQTL. For instance, the MQTL2B.3 overlapped with 132 MTAs, followed by the MQTL6B.1 which overlapped with 122 MTAs, while, several MQTLs each overlapped with only one MTA. Some MQTLs (e.g., MQTL1B.3, MQTL2B.3, MQTL2D.1, and MQTL4B.2) over-   Fig. 3). These genes encode different types of proteins involved in different pathways, for instance, transcription factors (e.g., AS2/LOB, ERF, MYB, NAC, and WRKY), protein kinases (e.g., CBL-interacting protein kinase, LRR receptor-like kinase protein, phosphoenolpyruvate carboxylase kinase), transporters (e.g., aluminiumactivated malate transporter) and other proteins including ABA 8´-hydroxylase, ABI5 binding protein, allene oxide cyclase, aminocyclopropane-1-carboxylate oxidase, GA-stimulated transcript family protein, glyoxalase, nicotianamine synthase, etc.
Candidate genes available from MAST-MQTLs. As many as 539 gene models were detected in the genomic regions of 28 selected MQTLs, each based on at least nine QTLs (providing tolerance to at least five different abiotic stresses) and with an average genetic and physical CI of 0.65 cM and 59.66 Mb, respectively (Table S7). These gene models included 66 gene models with no function descriptions. At the two extremes, a maximum of 45 gene models were available from MQTL1D.4, whereas, a minimum of only 4 gene models were   In silico expression analysis of genes and common cis-regulatory elements in the promoters of multiple abiotic stress responsive genes. In silico expression analysis was performed for all 539 gene models available from 28 selected MQTLs. The first transcriptomic dataset (which included data on genes expressed under DS, HS, and combined DS and HS conditions) identified 70 differentially expressed genes (DEGs), including 25 up-regulated genes, 7 down-regulated genes, and 38 genes that were up-regulated in some conditions but down-regulated in others. MQTL7D.1 had the most DEGs (7)  Overall, 189 DEGs were identified from 27 MQTLs (no DEG was available from MQTL2A.5). Furthermore, 27 genes were common between the first and second transcriptomic datasets, and 36 genes were common between the third and fourth transcriptomic datasets. There were 11 genes which were common among first, second and third datasets; these genes included 7 genes encoding for previously characterized proteins (viz., late embryogenesis abundant protein, sugar transporter-like, O-acyltransferase, zinc finger protein, transmembrane protein, fatty acid hydroxylase, and protein kinase), and the remaining 4 genes encoding for uncharacterized proteins (including protein of unknown function DUF1666, protein of unknown function DUF247, uncharacterized conserved protein UCP031279). Among these 11 genes, two genes (encoding for sugar transporter-like and transmembrane proteins) showed differential expression across all the four transcriptomic datasets investigated ( Table 2). A representative heat map of these 11 most promising DEGs is presented in Fig. 4.
Gene control is mediated by cis-regulatory elements (CREs) in promoter regions. It is crucial to study the CREs found in the putative promoter regions of the DEGs reported in this study. CREs were studied in 1.5 kb 5′ upstream promoter regions of each of the 11 multiple abiotic stress-responsive genes. A total of 359 CREs were discovered, including 19 unnamed/uncharacterized CREs (Table S8). The CREs were divided into three categories: plant growth and development, stress responsiveness, and phytohormone responsiveness, in addition to the conventional cis-elements found in promoters of DEGs (Fig. 5). Shoot and root meristem expression (CAT-box), endosperm expression (GCN4-motif and AAGAA-motif), flowering (CCAAT-box, AT-rich element, and Circadian), and zein metabolism (O2 site) were all included in the plant growth and development category. Drought-inducibility (MBS), anaerobic induction (ARE), stress (TC-rich repeats, As-1 element, STRE), low temperature (LTR), and wounding responsiveness (WRE3 and WUN-motif) were all included in the stress responsiveness group. Many cis-elements, including salicylic acid (TCA-element), ethylene (ERE), Me-JA (TGACG-motif and MYC), auxin (TGA-element), and abscisic acid (ABRE), were shown to be related to phytohormone responsiveness. ABRE, MYC, and ERE were the most common cis-elements, and they were all related to hormone responsiveness. These findings show that complex regulatory networks play a role in the transcriptional regulation of genes that play a role in MAST.
Conserved genomic regions associated with multiple abiotic stress tolerance in different cereals. The syntenic relationship was also observed for wheat genes available from MQTL regions with the maize and rice genomes. The 5161 and 5893 wheat genes available from MQTLs were used to detect the 5842 maize orthologues, and 5277 rice orthologues, respectively (Fig. 6a,b). This synteny analysis revealed several genomic regions in both maize and rice genomes where several MQTLs have already been identified in previous studies for different abiotic stresses; such genomic regions conserved between 'wheat and maize' and 'wheat and rice' may be termed as 'ortho-MQTLs' . Out of 5842 maize genes, 2565 genes were available in 101 maize MQTL regions previously identified to be associated with drought stress tolerance [40][41][42] . Similarly, out of 5277 rice genes, 2704 genes belonged to 140 rice MQTLs earlier reported to be associated with heat, drought, and salinity stress tolerance [43][44][45][46][47] . The number of conserved genes at corresponding orthologous positions of MQTLs (between 'wheat and maize' and 'wheat and rice') ranged from just a few to hundreds of genes. Detailed information related to maize and rice orthologues of wheat genes and ortho-MQTLs is presented in Table S9. Some earlier

Discussion
DS, HS, SS, WS, PHS, and AS are all abiotic stresses that are detrimental to wheat, posing significant intimidation to the environment and agriculture and resulting in considerable crop yield loss. These abiotic stresses may be avoided by growing wheat varieties conferring tolerance to individual or multiple abiotic stresses. Thus, new  www.nature.com/scientificreports/ wheat cultivars must contain novel loci/genes to provide tolerance against these stresses that influence plant development and yield. A number of studies on QTL mapping for different abiotic stresses in wheat have been undertaken over the last two decades (Table S1). Most of the QTLs found in these previously reported studies had a long CI and a relatively low PVE, which made them ineffective for use in MAB to develop wheat varieties providing tolerance against different abiotic stresses. Furthermore, validation is also required whenever these QTLs are planned to be used in a wheat breeding program. In comparison to the QTLs that have been previously identified, MQTLs are observed to be more stable and consistent, with a narrow CI and comparatively higher confidence, making them useful for crop improvement programs as well as for basic research programs, including the cloning and further characterization of favorable QTLs/genes 26 . Furthermore, several QTLs have been detected and characterized, each conferring tolerance to multiple abiotic stresses in different crops, including wheat 10,49-54 . Furthermore, numerous studies have uncovered different genes each conferring tolerance to multiple abiotic stresses in wheat, which are regarded as extremely valuable assets for future breeding programs aimed at developing climate-smart wheat varieties 14,15,55 . Based on these studies, the MAST hypothesis was figured out, which says that one locus is believed to make a plant more resistant to more than one abiotic stress. Several meta-QTL analyses for various traits in cereal crops such as maize 41 , rice 20,43,48 , and barley have already been conducted 35 . In wheat, meta-QTL studies have been undertaken for four different abiotic stresses, including PHS 56 , SS 30 , DS 32 , and HS 33 , but no meta-analysis has been conducted for the traits contributing to WS and AS stresses. In addition, the information from such studies quickly becomes outdated due to the identification and characterization of new QTLs in recently conducted studies. This necessitates a comprehensive meta-QTL analysis at a regular interval for proper utilization in the breeding pipeline. Furthermore, to our knowledge, no MAST meta-analysis has ever been conducted in wheat. In this study, we investigated the QTLs for six different abiotic stresses that have been reported until 2021.
This consensus map developed in the present study is much better and more comprehensive than previous consensus maps used for meta-analysis for some of these abiotic stress tolerance traits 30,32,33 . The consensus map utilized in the present study contained 100,614 markers (mostly SSR and SNPs) scattered throughout 6647 cM. Based on the QTL information and the development of the consensus map in this study, it is observed that the sub-genome B carries the maximum number of markers and contains the highest density of QTLs, while the sub-genome D has the minimum number of markers and contains the lowest density of QTLs. The characteristics of marker distribution and density of QTLs on different wheat chromosomes analyzed in the current meta-study are in agreement with other meta-analyses studies on individual abiotic stress tolerance earlier conducted in wheat [31][32][33] .
In this study, we predicted 134 MQTLs using previously identified QTLs. The majority of MQTLs (87.31%; 117/134) were each associated with tolerance against at least two of the six stresses under study. Furthermore, 23 MQTLs each were involved in providing tolerance against at least five different stresses, while five MQTLs, each were associated with all six abiotic stresses under study. As discussed, several MQTLs contained QTLs for different abiotic stress tolerance traits, suggesting a strong association between different abiotic stress tolerance traits. Furthermore, QTLs for the same abiotic stress tolerance from different mapping populations were projected onto the same chromosomal region, confirming that those regions exist. The meta-analysis conducted in this study refined the CIs and compared the genomic positions of different QTLs detected in individual studies associated Further, GWAS is an alternative high-precision technique to QTL mapping which is based on high-throughput sequencing or array technologies [57][58][59] . Interestingly, more than 57% (76/134) of MQTLs were validated with GWAS-MTAs associated with different abiotic stress tolerance traits identified in wheat. The MQTL2B.3 overlapped with the maximum number of MTAs (132) obtained from multiple GWA studies, while some MQTLs overlapped with only one MTA. Some recent studies also provided similar results, such as 54.6%, 63%, and 61.3% verification of MQTLs with GWAS, respectively 23,26,27 . Furthermore, the MQTLs validated using GWAS data are important genomic regions that could be utilized effectively in future wheat breeding programs and cloning of MAST genes available from MQTLs.
A comparison of MQTLs with known genes for distinct abiotic stresses in wheat could aid investigation into the molecular mechanisms that confer tolerance to various abiotic stresses. As a result, the relationships between MQTLs and known genes were also investigated in the present study. Twenty-one MQTLs were found to be co-localized with 34 genes in wheat that are involved in providing tolerance to different abiotic stresses. For instance, MQTL2B.6 overlapped with 2 PHS tolerance genes (TaSdr-A1 35 and TaAFP-B1 60 TaNAC67 17 ). According to these findings, the twenty-one MQTLs, which overlapped with 34 known genes for different abiotic stresses, appear to be promising candidates for improving tolerance to multiple abiotic stresses and could be valuable assets in breeding for climate resilience. However, more work is required to find functional variants of these genes in the MQTL regions that are being studied. Furthermore, there are several genes known in wheat that confer tolerance to multiple abiotic stresses 16,35,65,72,73 . Unfortunately, no MAST gene or QTL providing tolerance against all the six abiotic stresses under this study has been identified. In the present study, mining of genes in genomic regions of the 28 selected MAST-MQTLs allowed the identification of 539 non-redundant gene models. In silico expression analysis unraveled the 11 promising CGs showing differential expressions in at least three transcriptomic datasets. These genes included 7 genes encoding for characterized proteins and the remaining 4 genes encoding for uncharacterized proteins. These genes exhibited differential expressions in a temporal and stress-specific (or non-specific) manner when subjected to different abiotic stresses, for instance, gene TraesCS2B02G319600 exhibited significant up-regulation after one hour of DS treatment, whereas, it showed significant down-regulation after 6 h of DS treatment indicating that there is a rapid response to DS at the molecular level via translational regulation. On the other hand, genes TraesCS2D02G211100, TraesCS3B02G040300, and TraesCS3B02G466000 consistently showed up or down-regulation of transcript levels under HS and DS conditions, regardless of the time point analyzed; this indicates that these genes might be involved in multiple signal pathways in adaptation to DS and HS, as reported in some earlier studies 74,75 .
The relevance of some of these genes encoding for known proteins can be summarized as follows: (i) Late embryogenesis-abundant (LEA) proteins are a large and diverse family of proteins that are hypothesized to play a role in normal plant growth and development as well as protecting cells against a variety of abiotic stresses. Chen and colleagues 76 identified and characterized several LEA family genes conferring tolerance to multiple environmental stresses 76 . (ii) In a recent study, Xuan and colleagues 77 reported that a higher expression level of sugar transporter genes regulates plant responses to different stresses, including drought, salt, and low-temperature stress. (iii) Transgenic plants with increased expression of the WSD1 gene, which encodes for wax synthase/acyl-CoA:diacylglycerol acyltransferase showed improved tolerance to abscisic acid, mannitol, drought, and salinity 78 . These findings showed that manipulating cuticular waxes could be beneficial for increasing plant productivity in a rapidly changing climate. (iv) Overexpression of ZFP182, a TFIIIA-type zinc finger protein significantly improved tolerance against salt, cold, and drought stresses in transgenic rice 79 . Further, zinc finger proteins have been termed the master regulators of abiotic stress responses in plants 80 . (v) Over-expression of the A. thaliana SOS1 gene, which encodes a transmembrane protein, improves plant salt tolerance in A. thaliana 81 . (vi) In A. thaliana, AtBI-1 gene plays important roles in the resistance to various biotic and abiotic stresses. Studies have demonstrated that sphingolipid fatty acid 2-hydroxylase is required for the function of AtBI-1 gene 82 . However, these genes are assumed to have a pleiotropic function, but there is still no evidence of how these genes may confer tolerance to multiple abiotic stresses; further studies are required to characterize their functionality in response to multiple abiotic stresses. Moreover, the DEGs with no function descriptions available may be characterized for their possible roles in MAST in wheat in further studies.
Since there is a cross-talk between different abiotic stresses, identification of common and/or overlapping CREs is crucial for developing wheat cultivars that show tolerance towards them. CREs are DNA sequences found in gene promoters that transcription factors bind to, and they further drive gene expression by either enhancing or silencing transcription regulations in response to unfavorable environmental conditions 83 . In the present study, all 11 multiple abiotic stress-responsive DEGs were identified to be equipped with three different functional categories of cis-elements, such as plant growth and development, stress responsiveness, and phytohormone responsiveness. In the case of stress responsiveness cis-elements, for instance, the MYB binding site (MBS) element present in the ZmSOPro promoter region triggers drought stress tolerance in transgenic Arabidopsis 84 . Stress-responsive gene regulation in Arabidopsis and soybean essentially depends on the TGG TTT (ARE), activation sequence-1 (as-1) and TC-rich repeat elements [85][86][87] . In addition, WUN-motif is a dehydration-related element detected in the regulatory region of the StDREB gene in tomato 88 .
The current study additionally investigated chromosomal regions with MQTLs that are conserved across several cereals. These ortho-MQTL regions were identified by integrating the two different approaches utilized by Saini et al. 26 and Sandhu et al. 20 , respectively. Ortho-MQTL analyses have recently been conducted in wheat for different traits including quality traits 28 , HS 33 , SS 30 , nitrogen use efficiency 24 , grain yield and component traits 22 . The discovery of ortho-MQTLs in closely related species shows their relevance and stability, as well as the reliability of associated genes 56 . Functional validation of the orthologous genes uncovered in this study is feasible, and gene-based functional markers for cereal breeding programs can be developed 89 .
Further, based on the criteria given in Saini et al. 26  www.nature.com/scientificreports/ and environments. They can be used either in MAB programs aimed at developing wheat varieties providing tolerance to multiple abiotic stresses under study or subjected to fine mapping to identify key genes involved in tolerance against multiple abiotic stresses in wheat, as proposed recently for multiple biotic stress resistance in soybean 101 . Markers flanking these MQTLs may be utilized in improving the prediction accuracy of genomic selection for multiple abiotic stress tolerance in wheat 102,103 .

Concluding remarks
The results of interval mapping studies on DS, HS, SS, WS, PHS, and AS tolerance were successfully integrated during the current study, resulting in the identification of 134 MQTLs. GWAS was used to confirm more than half of these MQTLs. Thirty-four known genes associated with various abiotic stresses overlapped several of these MQTLs as well. Furthermore, 28 MQTLs conferred tolerance to five or all of the six abiotic stresses studied. Putative genes associated with several abiotic stresses were discovered, and 11 genes encoding important proteins showed differential expression across different transcriptome datasets. For the development of tolerance to numerous abiotic stresses in wheat cultivars, five promising MQTLs imparting tolerance to at least five abiotic stresses were recommended for use in marker-assisted breeding. Overall, these findings provide useful information as well as robust candidate genes for further functional investigation aimed at enhancing wheat MAST.

Material and methods
Collection of QTL mapping studies and data on QTLs.  104 ). Each independent study was further used to obtain the following information: (i) markers flanking the QTLs, (ii) information on mapping populations (size and type), (iv) genetic peak position and confidence interval (CI), (v) LOD score, and (vi) phenotypic variation explained by individual QTLs (R 2 or PVE). Studies with incomplete datasets were excluded. In some cases, the peak position of the QTLs were not available, consequently, the mid-point of both the flanking markers were considered the peaks. Simultaneously, a LOD score of 3.0 and a PVE value of 10 were assigned wherever no LOD score or PVE value was provided. All of the QTLs involved in this study were given unique identities as follows: the author's name followed by the year of publication, the abbreviated name of the trait, stress, and chromosome involved, respectively. Numerals after the chromosome numbers were used to distinguish various QTLs for the same trait and stress on the same chromosome identified in the same study. Under each stress condition, QTLs associated with different traits, including yield and related traits, and several morpho-physiological traits, were reported in the literature. For the purpose of this analysis, all the QTLs associated with different component traits identified under a particular stress condition were treated as QTLs associated with the abiotic stress in question.

Construction of the consensus map and projection of QTL.
A highly dense consensus map was constructed by integrating the following high-quality linkage maps: (a) 'Wheat, Consensus SSR, 2004′ 105 (b) 'Synthetic x Opata BARC map' 106 and SNP maps assayed through (c) 'Illumina iSelect 90 K SNP Array' 107 , (d) 'Wheat 55 K SNP array' 108 . In addition, markers flanking the individual QTLs were also integrated into this consensus map to ensure the projection of the maximum number of QTLs. The R package 'LPMerge' was utilized to generate the consensus map 109 .
Separate map and QTL data input files from several QTL mapping experiments were prepared. The map file included information mainly on linkage groups, markers, and their genetic positions, whereas, the QTL file contained information on QTL name, stress name, chromosome number, LOD, PVE, and genetic positions (peak position and CI) of individual QTLs. In some studies where the original CI for individual QTLs was not reported, CI (95%) was obtained for each such QTL using different population-specific equations [110][111][112] . According to Chardon et al. 113 QTLs were projected onto the newly created consensus map using the projection tool (QTLProj) in BioMercator v4.2.3 software 18 (https:// sourc esup. renat er. fr/ forum/ forum. php? forum_ id= 5025). To select the optimum projection context, QTLProj employs a dynamic strategy. A perfect context consists of a pair of common markers that flank the QTL in the original map and a consistent distance estimate between the two maps (initial and consensus). The behaviour of the procedure to get such a configuration is influenced by the minimal value of the flanking marker interval distances ratio and the minimal p-value produced by confirming the homogeneity of the flanking marker interval distances between the initial and consensus maps 19 . Prediction of MQTLs and their validation using GWAS. Individual MQTLs were predicted for each wheat chromosome as per instructions provided in the manual (https:// www. ebi. ac. uk/ eccb/ 2014/ eccb14. loria. fr/ progr am/ id_ track/ ID10-summa ry. pdf) and keeping abiotic stress as a meta-trait (including all six abiotic stresses: DS, HS, SS, WS, PHS, and AS). The two-steps method proposed by Veyrieraset al. 19 was utilized for the prediction of MQTLs. The first step involves the selection of the best MQTL model based on obtaining the lowest values of the selection criteria in at least three of the following five models: (i) AIC (Akaike information content), (ii) AICc (AIC correction), (iii) AIC3 (AIC 3 candidate models), (iv) BIC (Bayesian information criterion), and (v) AWE (average weight of evidence). The second step predicts the total number of MQTLs on a chromosome along with information about peak position, CI (95%), and weightage based on the model selected. An MQTL was described as a genomic region in which two or more QTLs from distinct mapping studies were involved. The means of the LOD and PVE values of the initial QTLs involved were used to compute the LOD and PVE values of the predicted MQTLs. When more than one MQTL was found on a single chromosome, the MQTLs www.nature.com/scientificreports/ were named using the standard procedure, which begins with the individual identification of each chromosome, followed by a point and a number (e.g., MQTL1A.1, MQTL1A.2). The physical positions were obtained for each MQTL using the available sequences of each marker flanking the MQTLs. These nucleotide sequences were retrieved for each marker, using either the published studies or online databases, including GrainGenes (https:// wheat. pw. usda. gov/ GG3/), Gramene (https:// www. grame ne. org/), and CerealsDB (https:// www. cerea lsdb. uk. net/ cerea lgeno mics/ Cerea lsDB /indexNEW.php). Consequently, these nucleotide sequences were subjected to BLASTN searches against the reference genome of wheat available in the EnsemblPlants database (https:// plants. ensem bl. org/ index. html). Genomic regions of some SNPs were obtained from the URGI-JBrowse wheat genome browser (https:// wheat-urgi. versa illes. inra. fr/). Furthermore, data from 32 independent GWA studies (involving 10 GWA studies on DS, 4 on HS, 9 on SS, 2 on PHS, and 1 on AS, 6 studies including both DS and HS; no GWAS study was available for WS) were collected to validate the MQTLs discovered in the current study (Table S1). GWA studies, involving durum wheat, hexaploid wheat (e.g., spring and winter wheat), and synthetic hexaploid wheat were conducted across 11 countries. The population size ranged from 70 114 to 717 115 . The population size, abiotic stress associated, genotyping platform, number of SNP markers employed, and significant markers identified in individual studies are all listed in Table S1. The physical positions of identified significant markers were retrieved from source studies or databases (such as the URGI-JBrowse wheat genome browser and CerealsDB) and compared to the physical coordinates of MQTLs identified on individual chromosomes. MQTLs co-localizing with at least one significant marker were accepted as GWAS validated/verified MQTLs.
Co-localization of known genes with MQTLs. The association between MQTLs and known genes involved in stress tolerance was also investigated. Sequences of these genes or related markers were BLASTed against the wheat reference genome (IWGSC RefSeqv1.0) in the EnsemblPlants database to determine their physical positions. The physical positions of known genes were compared with the physical positions of MQTLs to identify the co-localization of MQTLs with the known genes.
Candidate gene mining, expression analysis, and promoter analysis. MQTLs each conferring tolerance to five or all six abiotic stresses were selected for the identification of underlying gene models. Of the selected MQTLs, MQTLs with a physical CI of ≤ 2 Mb were directly utilized for gene mining. Only a 2 Mb physical region (1 Mb physical region on either side of the MQTL peak) was examined for the availability of gene models for the remaining MQTLs with longer physical CI. The peak genomic positions of the MQTLs were estimated using the formula proposed by our group recently 26 , which is as follows: The InterPro database (https:// www. ebi. ac. uk/ inter pro/) was used to obtain descriptions of the functions of discovered gene models. In silico expression analysis of the available gene models was performed using transcriptome data from the wheat expression browser (expression and visualisation platform) powered by expVIP (http:// www. wheat-expre ssion. com116). For this purpose, the following relevant expression datasets were utilized: (i) 'drought and heat stress time-course in seedlings' 74 , (ii) 'seedlings with PEG to simulate drought' 116 , (iii) 'spikes with water stress' 116 , and (iv) 'grain developmental time-course with 4A dormancy QTL' 117 . The first dataset 74 included differential gene expression in the drought and heat tolerant wheat cultivar TAM107 grown under DS (20% PEG-6000), HS (400C), and both DS and HS conditions (40 °C and 20%, PEG-6000), with leaf samples collected separately at 1 and 6 h after treatment. Seedlings in normal growth conditions (well-watered, 22 °C) were used as controls in all of the experiments. The second dataset 116 consisted of differential expression of genes in two cultivars (Gemmiza10 and Giza168) grown under DS (PEG 6000) with samples from seedlings collected 2 and 12 h after treatment. The third dataset 114 consisted of differential expression of genes in a DH population (Westonia x Kauz) grown under water stress conditions with samples collected from the spike at the early booting stage. The fourth dataset 117 consisted of expression data from 9 NILs (segregating for the 4A dormancy QTL region) and four parental varieties (viz., Baxter, Chara, Westonia, and Yipti) with grain samples collected as 15, 25, and 35 days post-anthesis. The remaining three abiotic stresses, SS, WS, and AL, did not have expression datasets accessible from expVIP. The expression data was given as log2 transformed TPM (transcripts per million) values. Only genes with a substantial fold change (FC ≥ 2 or FC ≤ − 2 TPM values) compared to control were considered differentially expressed. To represent gene expression patterns, heat maps were created using the web tool Morpheus (https:// softw are. broad insti tute. org/ morph eus/). Further, the 1500-bp nucleotide sequences upstream of the translational start site (i.e., ATG) were retrieved from the EnsemblPlants database and uploaded to the PlantCARE database 118 for cis-regulatory elements (CREs) prediction and potential activities. Only response elements available only on the sense strand with a matrix value of 5 were acceptable, as reported in Sharma et al. 119 .
Revealing conserved genomic regions associated with multiple abiotic stress tolerances in different cereals. Some of the most promising MQTLs (each conferring tolerance to 5 or all 6 stresses in question) were selected for identification of conserved genomic regions (including ortho-or syntenic MQTL regions) among the wheat, rice, and maize genomes. This analysis involved the following two steps: (i) gene models from the selected MQTL regions were BLASTed against rice and maize genomes available at EnsemblPlants database to discover conserved genomic regions, and (ii) the matching rice and maize orthologues were retrieved from the database with their genomic coordinates. In order to identify the ortho-MQTLs from the genomic